
clear all;
close all;

MAX_NUM_ITER = 50000;

A = cell(4, 1);

# Newton-based algorithms testing the Wolfe conditions
A{1}.algorithms = { "gsl_ext_bfgs_mt", "gsl_ext_bfgs_f", "gsl_vector_bfgs2", "gsl_ext_mnewton" };
A{1}.name = "newtonbased_wolfe";

# conjugate gradient algorithms testing the Wolfe conditions
A{2}.algorithms = { "gsl_ext_conjgrad_fr_mt", "gsl_ext_conjgrad_pr_mt" };
A{2}.name = "conjrad_wolfe";

# gradient descent algorithms with the Brent line search
A{3}.algorithms = { "gsl_vector_bfgs", "gsl_conjugate_fr", "gsl_conjugate_pr" };
A{3}.name = "gsl_brent";

# simplex algorithms
A{4}.algorithms = { "gsl_ext_lrwwsimplex", "gsl_nmsimplex" };
A{4}.name = "simplex";

testfunctions = GSLpp_testfunctions();

for tfi = 1:length(testfunctions)
  if strcmp(testfunctions{tfi}.stopcrit_name, "xdisttomin") == 0
    continue;
  endif

  for asi = 1:length(A)
    R = NaN(MAX_NUM_ITER, length(A{asi}.algorithms));
    maxnumiter = 0;
    maxnorm = 1e6;

    for ai = 1:length(A{asi}.algorithms)
      tf_params = GSLpp_get_testfunction_params(testfunctions{tfi});

      if strcmp(typeinfo(testfunctions{tfi}.starting_point), "function handle") == 1
        x0 = feval(testfunctions{tfi}.starting_point, testfunctions{tfi}.n);
      else
        x0 = testfunctions{tfi}.starting_point;
      endif

      stopcrit_params = testfunctions{tfi}.stopcrit_params;
      stopcrit_params.eps = 1e-12;

      results = GSLpp_minimize(A{asi}.algorithms{ai}, ...
                struct, struct, ...
                testfunctions{tfi}.name, tf_params, ...
                struct("f", "sym", "g", "sym", "H", "sym"), ...
                "xdisttomin", stopcrit_params, ...
                x0, 0, false);

      AR = GSLpp_convergence_plot(results, { 'disttominimizer' });
      AR = reshape(AR, length(AR), 1);

      if results.numiter < MAX_NUM_ITER
        maxnumiter = max(maxnumiter, results.numiter);
      endif
      maxnorm = min(maxnorm, results.iterations{1}.stopcrit_value);

      R(1:results.numiter, ai) = AR;
    endfor

    if maxnumiter == 0
      maxnumiter = MAX_NUM_ITER;
    endif

    R = resize(R, maxnumiter * 1.2, size(R, 2));

    L = A{asi}.algorithms;
    for j = 1:length(L)
      L{j} = GSLpp_title(L{j});
    endfor

    GSLpp_newplot();
    semilogy(R, "linewidth", 3);
    legend(L);
    legend("boxon");
    title(testfunctions{tfi}.title);
    xlabel("Iterations");
    ylabel("||x_{k}-x^{*}||");
    axis([1, length(R), 1e-12, maxnorm]);

    print([ "convrates_", A{asi}.name, "_", testfunctions{tfi}.name, ".ps" ], "-landscape", "-dashed", "-FArial:23");
    close;
  endfor
endfor
